import numpy as np
import cv2
import matplotlib.pyplot as plt

# 读取图像并手动转换为RGB格式
img = cv2.imread('mytest1.jpg')
img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)


# Sobel算子滤波
def sobel_filter(img):
    height, width = img.shape[:2]
    result = np.zeros_like(img)

    # Sobel算子
    sobel_x = np.array([[-1, 0, 1],
                        [-2, 0, 2],
                        [-1, 0, 1]])
    sobel_y = np.array([[-1, -2, -1],
                        [0, 0, 0],
                        [1, 2, 1]])

    # 对每个通道进行卷积
    for c in range(3):
        # 手动填充边缘
        pad_img = np.zeros((height + 2, width + 2))
        pad_img[1:-1, 1:-1] = img[:, :, c]
        pad_img[0, 1:-1] = img[0, :, c]
        pad_img[-1, 1:-1] = img[-1, :, c]
        pad_img[1:-1, 0] = img[:, 0, c]
        pad_img[1:-1, -1] = img[:, -1, c]
        pad_img[0, 0] = img[0, 0, c]
        pad_img[0, -1] = img[0, -1, c]
        pad_img[-1, 0] = img[-1, 0, c]
        pad_img[-1, -1] = img[-1, -1, c]

        for i in range(height):
            for j in range(width):
                # x方向梯度
                gx = 0
                for k in range(3):
                    for l in range(3):
                        gx += pad_img[i + k, j + l] * sobel_x[k, l]
                # y方向梯度
                gy = 0
                for k in range(3):
                    for l in range(3):
                        gy += pad_img[i + k, j + l] * sobel_y[k, l]
                # 梯度幅值
                result[i, j, c] = np.sqrt(gx ** 2 + gy ** 2)

    return np.uint8(result)


# 给定卷积核滤波
def custom_filter(img):
    height, width = img.shape[:2]
    result = np.zeros_like(img)

    # 给定卷积核
    kernel = np.array([[1, 0, -1],
                       [2, 0, -2],
                       [1, 0, -1]])

    # 对每个通道进行卷积
    for c in range(3):
        # 手动填充边缘
        pad_img = np.zeros((height + 2, width + 2))
        pad_img[1:-1, 1:-1] = img[:, :, c]
        pad_img[0, 1:-1] = img[0, :, c]
        pad_img[-1, 1:-1] = img[-1, :, c]
        pad_img[1:-1, 0] = img[:, 0, c]
        pad_img[1:-1, -1] = img[:, -1, c]
        pad_img[0, 0] = img[0, 0, c]
        pad_img[0, -1] = img[0, -1, c]
        pad_img[-1, 0] = img[-1, 0, c]
        pad_img[-1, -1] = img[-1, -1, c]

        for i in range(height):
            for j in range(width):
                sum = 0
                for k in range(3):
                    for l in range(3):
                        sum += pad_img[i + k, j + l] * kernel[k, l]
                result[i, j, c] = sum

    return np.uint8(result)


# 计算颜色直方图
def color_histogram(img):
    hist = np.zeros((3, 256))
    height, width = img.shape[:2]

    for c in range(3):
        for i in range(height):
            for j in range(width):
                hist[c, img[i, j, c]] += 1

    # 手动归一化
    total = height * width
    for c in range(3):
        for i in range(256):
            hist[c, i] = hist[c, i] / total

    return hist


# 提取纹理特征 (使用灰度共生矩阵)
def texture_features(img):
    # 手动转换为灰度图
    height, width = img.shape[:2]
    gray = np.zeros((height, width), dtype=np.uint8)
    for i in range(height):
        for j in range(width):
            gray[i, j] = 0.299 * img[i, j, 0] + 0.587 * img[i, j, 1] + 0.114 * img[i, j, 2]

    glcm = np.zeros((256, 256))

    # 计算水平方向的灰度共生矩阵
    for i in range(height):
        for j in range(width - 1):
            glcm[gray[i, j], gray[i, j + 1]] += 1

    # 手动归一化
    total = np.sum(glcm)
    for i in range(256):
        for j in range(256):
            glcm[i, j] = glcm[i, j] / total

    # 计算特征
    contrast = energy = correlation = homogeneity = 0

    # 手动计算均值和标准差
    mean = 0
    for i in range(height):
        for j in range(width):
            mean += gray[i, j]
    mean = mean / (height * width)

    std = 0
    for i in range(height):
        for j in range(width):
            std += (gray[i, j] - mean) ** 2
    std = np.sqrt(std / (height * width))

    for i in range(256):
        for j in range(256):
            contrast += (i - j) ** 2 * glcm[i, j]
            energy += glcm[i, j] ** 2
            correlation += (i - mean) * (j - mean) * glcm[i, j] / (std ** 2)
            homogeneity += glcm[i, j] / (1 + abs(i - j))

    features = np.array([contrast, energy, correlation, homogeneity])
    return features


# 处理图像
sobel_result = sobel_filter(img)
custom_result = custom_filter(img)
hist = color_histogram(img)
texture = texture_features(img)

# 保存结果
cv2.imwrite('sobel_result.jpg', cv2.cvtColor(sobel_result, cv2.COLOR_RGB2BGR))
cv2.imwrite('custom_result.jpg', cv2.cvtColor(custom_result, cv2.COLOR_RGB2BGR))
np.save('texture_features.npy', texture)

# 显示结果
plt.figure(figsize=(15, 10))

plt.subplot(231)
plt.imshow(img)
plt.title('原图')

plt.subplot(232)
plt.imshow(sobel_result)
plt.title('Sobel算子滤波结果')

plt.subplot(233)
plt.imshow(custom_result)
plt.title('给定卷积核滤波结果')

plt.subplot(234)
plt.plot(hist[0], 'r')
plt.plot(hist[1], 'g')
plt.plot(hist[2], 'b')
plt.title('颜色直方图')

plt.tight_layout()
plt.show()
